function [x, fv, error, steps] = min_Newton(func, dfunc, d2func, x0, MAX_IT, tol)
    err = 1.0;
    steps = 0;
    convergent = 1; 
    alpha = 1.0;
    x = [0, 0]';
    while (err > tol)
        f = func(x0(1), x0(2));
        g = dfunc(x0(1), x0(2));
        G = d2func(x0(1), x0(2));
        p = -G\g;   % Gp = -g
        x0 = x0 + alpha * p;
        steps = steps + 1;
        x(:,steps) = x0;
        fv(steps) = f;
        if (steps > MAX_IT)
            convergent = 0;
            break;
        end
        err = norm(g);
        error(steps) = err;
    end
    if (convergent == 0)
        steps = -1;
    end
end 